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1  Background 

The  preliminary  report  on  this  contract  (submitted  in  January  2002)  de¬ 
scribed  the  potential  of  microcavities  stabilizing  hypersonic  flows,  especially 
with  respect  to  the  second  mode  of  instability  -  the  dominant  source  of  flow 
transition  in  this  regime.  This  work  was  presented  at  the  3rd  AIAA  The¬ 
oretical  Fluid  Mechanics  Meeting  in  St  Louis  in  June,  and  an  offprint  of 
the  paper  (AIAA  2002-2987)  is  attached  as  an  appendix  to  this  final  report. 
The  first  phase  of  the  work  was  partly  concerned  with  the  nature  of  the 
flow  inside  each  individual  microcavity,  and  it  was  shown  that  this  could 
be  regarded  as  incompressible  and  highly  viscous.  This  result  was  exploited 
in  the  the  associated  work  of  Theofilis  (performed  under  EOARD/USAF 
Contract  No.  F61775-01-WE049)  which  focused  mainly  on  the  flow  inside  a 
microcavity  in  which  the  nature  of  the  basic  flow  was  modeled  as  that  inside 
a  single  cavity,  driven  by  a  uniform  shear  on  the  (upper)  boundary.  As  such, 
this  model  may  be  regarded  as  a  modification  of  the  classical  lid-driven  cav¬ 
ity  problem  of  Burggraf  [1],  with  a  single  modified  boundary  condition.  This 
work  also  indicated  that  the  flow  inside  each  microcavity  was  exceedingly 
stable,  i.e.  flow  instabilities  inside  each  microcavity  are  not  an  issue. 

Other  work  performed  in  the  first  phase  of  this  contract  focused  on  the 
(inviscid)  stability  of  the  outer  (boundary-layer)  flow,  modeling  the  wall 
porosity  by  means  of  a  modified  wall  boundary  condition,  linking  the  normal 
velocity  component  on  the  wall  to  the  wall  pressure.  This  indicated  strongly, 
that  there  certainly  exists  a  significant  potential  for  flow  stabilization,  over 
a  very  broad  range  of  parameters.  Typically  reductions  in  growth  rates 
of  circa  50%  seem  to  be  (easily)  possible,  and  this  is  in  surprisingly  good 
quantitative  agreement  with  the  experimental  work  of  Rasheed  et  al.  [2] 
(whilst  also  confirming  the  theoretical  work  of  [3] ) . 

^This  material  is  based  upon  work  supported  by  the  European  Office  of  Aerospace 
Research  and  Development,  Air  Force  Office  of  Scientific  Research,  Air  Force  Research 
Laboratory,  under  Contract  No.  F61775-01-WE046  monitored  by  Dr.  Steven  Walker 
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The  key  issue  addressed  in  the  second  phase  of  this  work  has  been  the 
linkage  between  the  flow  inside  the  microcavities  and  that  in  the  outer 
(boundary- layer)  flow;  this  final  report  details  results  related  to  this  aspect, 
and  also  indicates  the  direction  of  other  possible  studies  for  the  future. 


2  Description  of  work  performed 


The  outer  flow 


As  before,  taking  the  x  coordinate  in  the  streamwise  direction,  Y  the  normal 
boundary-layer  coordinate  (scaled  with  R  2  in  the  usual  boundary-layer 
manner),  then  the  corresponding  base  flow  (nondimensionalized  with  respect 
to  the  freestream  velocity  U0 c,  which  is  also  used  in  our  definition  of  the  outer 
Reynolds  number  R).  and  temperature  field,  then  takes  the  form 

(u,  T)  =  ((E/o(Y),0),T„(Y))  +  O(iH),  (1) 

where  the  use  of  the  parallel-flow  approximation  is  explicitly  assumed,  con¬ 
sistent  with  the  scale  of  the  flow  disturbance  field  (note  too  that  the  simi¬ 
larity  form  for  the  compressible  boundary-layer  equations  was  used  in  the 
evaluation  of  the  above  quantities) . 

Small  amplitude  disturbances  are  imposed  on  (1)  and  written  in  the 
usual  (normal-mode)  form,  i.e. 

{u,v,p,T,p)  =  {u,av,p,f,p)etR*a(x-ct\  (2) 


Here  a  denotes  the  wavenumber  of  the  disturbance  and  c  is  the  associated 
wavespeed.  Substituting  the  disturbance  flowfield  into  the  inviscid,  com¬ 
pressible  equations  of  motion,  discarding  quadratic  amplitude  terms  leads 
to  the  following  system 
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The  above  may  be  reduced  to  two  first-order  equations,  namely 


Uqyv  =  JP_  \T0  -  Ml(U0  -  c)2' 
U0  -  c  jMl  [  U0-c  \  ’ 
ia2(Uo  —  c)v  _  py 

%  =~Wl' 


(8) 

(9) 


The  usual  zero  permeability  on  the  wall  condition  (v(Y  =  0)  =  0)  is  re¬ 
placed  by  the  following  condition  linking  perturbation  normal  velocity  to 
perturbation  pressure,  i.e. 


v  —  Ap  =  0.  (10) 

Focusing  on  the  particular  case  of  a  freestream  Mach  number  M <*,  =  6  the 
first  phase  of  this  work  (see  also  the  appended  AIAA  paper)  extensively 
investigated  the  effects  of  altering  A  (which  in  general  will  be  complex).  It 
was  found  that  values  of  A  lying  within  the  third  quadrant  were  particularly 
effective  at  reducing  maximum  growth  rates,  especially  that  of  the  second 
mode  of  instability. 

The  central  issue  is  then  how  to  link  these  outer-flow  features  to  the  flow 
inside  each  microcavity.  This  is  addressed  below. 

Disturbance  flowfield  inside  each  microcavity 

As  noted  above,  the  early  part  of  the  study  revealed  the  flow  inside  each 
microcavity  to  be  fundamentally  very  slow  (i.e.  incompressible  and  viscous). 
This  second  phase  of  the  work  has  been  involved  with  developing  a  model 
to  obtain  some  connection  between  the  inner  and  outer  flow  regimes.  The 
model  that  has  been  developed  in  this  second  phase  builds  upon  the  work  of 
Duck  [4]  on  oscillatory  lid-driven  cavities.  The  fundamentals  of  the  model 
used  inside  each  microcavity  are  as  follows:- 

•  The  flow  is  viscous  and  incompressible. 

•  The  flow  is  two-dimensional. 


•  The  usual  permeability  and  no-slip  conditions  are  imposed  on  the  bot¬ 
tom  and  two  sides  of  the  cavity  (which  is  assumed  to  be  rectangular 
in  shape). 

•  The  flow  is  driven  solely  by  the  disturbance  (time-harmonic)  flowfield, 
i.e.  the  outer  base-flow  boundary  layer  has  no  direct  effect  on  the  mi¬ 
crocavity  flow  (this  approximation  was  suggested  by  the  investigations 
in  the  first  phase  of  the  study). 
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The  forcing  (through  the  time-harmonic  component  of  the  flow)  com¬ 
prises  two  components: 


—  through  forcing  on  the  normal  velocity  on  the  upper  boundary; 
generally  this  in  the  form: 

v(x,  y  =  H,t)  =  ~vt0p  sm(2Trx)e^int  +  c.c.,  (11) 


—  through  forcing  on  the  tangential  velocity  on  the  upper  boundary; 
generally  this  was  imposed  in  one  of  two  forms 

Model  (a): 

u(x,  y  =  H,t)  =  ^utope^nt+<p)  +  c.c.,  (12) 

Model  (b): 

u(x,y  =  H,t)  =  -utopel^nt+(^  sin(7nr)  +  c.c..  (13) 

Here  Vtop  and  Utop  are  both  constants  (the  latter  implying  a  uniform 
tangential  velocity  across  the  top  of  the  cavity  in  the  case  of  model 
(a),  namely  (12)).  Note  that  here  is  directly  related  to  ac  in  the 
analysis  of  the  outer  flow. 


It  should  be  noted  that  (11)  was  chosen  specifically  such  that  v  =  0  at 
x  =  0  and  x  =  1  in  line  with  the  impermeability  constraint  away  from  the 
pores/microcavities,  and  note  too  the  zero  total  mass  flux  into  each  cavity. 
(12)  was  chosen  because  of  its  simplicity,  whilst  (13)  was  chosen  in  order  to 
satisfy  the  no-slip  condition  away  from  the  pores.  Defining  a  streamfunction 
ip,  such  that  (u,v)  =  (tpy,  —ipx)->  then  the  problem  (formally)  is  therefore  as 
follows: 


-ittV2xpt 


d{-ip,V2-ip) 

d{x,y) 


(14) 


In  the  above,  non-dimensionalisation  of  lengthscales  has  been  performed 
using  the  width  of  the  cavity,  and  as  a  consequence  0  <  x  <  1,  and  so 
the  vertical  scale  is  0  <  y  <  H.  Velocities  are  non-dimensionalised  with 
respect  to  the  normal  velocity,  vtop  and  this  is  also  used  in  the  definition  of 
the  microcavity  Reynolds  number,  Re  (along  with  the  cavity  width);  time 
is  non-dimensionalised  with  respect  to  the  cavity  width  and  vtop-  It  should 
be  emphasised  that  in  the  experiments  of  Rasheed  et  al.  [2]  (the  regime  of 
interest),  both  the  pore  spacing  and  width  were  considerably  smaller  than 
the  outer  boundary-layer  thickness. 
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Since  numerically  it  was  found  in  the  first  phase  of  the  study  that  9{ac} 
was  typically  numerically  small,  as  a  first  step  in  modeling  the  cavity  flow 
we  set  flj  =  0  (although  of  course  from  the  point  of  view  of  the  stability  of 
the  outer  flow,  the  finiteness  and  sign  of  9{ac}  is  crucially  important;  note 
that  |9{ac}|  <<  1  is  a  feature  common  to  many  flow  stability  studies). 

This  being  the  case,  the  flow  was  then  decomposed  harmonically  as  fol¬ 
lows  (following  Duck  [4]) 

OO 

1l’n(x,y)emt,  (15) 

n—— oo 

where,  if  ff  is  taken  to  be  real,  then  'ip  n  =  complex  conjugate{i/ir, }.  It 
was  also  found  useful  to  introduced  the  vorticity  components,  (n  defined  as 
follows: 

Cn  =  V2W  (16) 

Then  (16)  was  used  in  conjunction  with 
—iQn(n  — 

m 

The  appropriate  boundary  conditions  are 


( d'tprn  SCn— m  dlpm  d(,n— m  \  1  _2  * 

^  \  dx  dy  dy  dx  )  Re  n 


ipn{ o  <  X  <  l:y  =  0)  =  ipny( 0  <  X  <  l,y  =  0)  =  ipn(x  =  0,0  <  y  <  H)  = 
i’nxix  =  0,0  <  y  <  H)  =  ipn(x  =  1,0  <  y  <  H)  =  ipnx(x  =  1,0  <  y  <  H)  =  0, 

(18) 

whilst  on  the  top  (y  =  H ):  if  n  ^  ±1: 

ipn{0  <  x  <  l,y  =  H)  =  ipny( 0  <  x  <  l,y  =  H)  =  0,  (19) 

but  for  n  =  ±1  then: 

i>±i(0  <  x  <l,y  =  H)  =  ^  (cos(2?ra:)  -  1)  ,  (20) 

and  in  the  case  of  model  (a): 

^(0  <x<l,y  =  H)  =  (21) 

dy  2  vtop 

whilst  in  the  case  of  model  (b): 

(0  <  x  <  1,  y  —  H)  —  si h{ttx).  (22) 

oy  Z  Vt0p 
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The  system  (15)  was  truncated  (at  n  =  ±JV,  where  typically  N  =  4  was 
adequate  for  acceptable  accuracy)  and  then  solved  using  a  (second-order) 
finite-difference  approach;  although  these  are  well  known  for  losing  stability 
at  higher  Reynolds  numbers  (Re),  this  is  not  an  issue  in  the  parameter 
regimes  of  interest,  and  equally  the  accuracy  in  such  regimes  is  perfectly 
adequate,  on  relatively  modest  grids. 

The  final  (and  key)  step  is  then  the  determination  of  A  (or  its  equivalent) 
in  (10).  In  the  previous  study  of  Fedorov  et  al  [3]  (in  which  the  flow  details 
inside  each  microcavity  were  not  considered  in  the  same  detail  as  the  present 
study,  with  each  cavity  being  assumed  to  be  of  infinite  depth),  an  averag¬ 
ing  process  was  implemented.  Here  we  follow  a  similar  idea,  although  the 
particular  details  are  somewhat  different.  For  the  purposes  of  determining 
the  coupling  between  the  wall  normal  velocity  and  the  pressure,  we  focus 
on  just  the  primary  (n  =  ±1)  modes.  The  system  (16)  together  with  (17) 
admits  the  following  symmetries: 


P2n+l(^  ~X,y)  =  -P2n+l(x~  ^,y),  1p2n+l  (^  ~  X,y)  =  tWl  (*  “  V) , 

j  i  .  ,  (23) 

P2n(^~  X,y)  =  p2n(x  ~  ~,y),  $2n(-  -  X,  y)  =  ~lfan(x  ~  - ,y).  (24) 

Here  pn  denotes  the  nth  pressure  mode.  First  an  averaging  process  was 
performed  across  the  top  of  each  cavity,  to  yield 

rl  ~fpix(x,y  =  H) 


A  = 


/ 


Pi(x,y  =  H) 


-dx. 


(25) 


In  order  to  evaluate  p\,  first  p\x  was  determined  from  the  x-component  of 
the  Navier-Stokes  equations: 


1 


pix  —  niipiy  + 


d^rn  d2^n- 
dx  dy'2 


m=— oo 

The  antisymmetry  (23)  then  enables  us  to  write 

rx  Qpx 


dlpm  d2V>n-?7 
dy  dydx 


Pi(x,y)  = 


dx 


dx. 


(26) 


(27) 


Finally  the  following  was  taken  as  an  estimate  for  the  coupling  coefficient 


A*  =  crA, 


(28) 


where  a  is  a  measure  of  the  porosity  of  the  surface,  which  is  the  fraction 
of  the  overall  volume  taken  up  by  the  pores;  following  Fedorov  et  al  [3]  we 
generally  took  a  =  0.5. 
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Vtop 


Figure  1:  Singular  case,  fl  =  1,  Re  =  20,  </>  =  7r/2  (square  cavity) 


Results  -  model  (a) 

A  number  of  calculations  have  been  performed  using  this  model.  It 
should  be  noted  that  in  this  case  a  singularity  on  the  top  wall  is  provoked  at 
x  =  0  and  x  =  1,  due  to  the  discontinuity  in  x-wise  velocity,  in  line  with  an 
identical  form  of  singularity  encountered  in  the  classical  steady  lid-driven 
case  of  Burggraf  [1]. 

A  number  of  preliminary  calculations  were  performed  with  various  values 
of  the  phase  shift  parameter  (j).  with  a  view  to  determining  the  possibility  of 
achieving  values  of  the  transpiration  parameter  A*  within  the  third  quadrant 
-  the  zone  found  in  the  first  phase  of  this  work  that  produced  the  best 
reductions  in  growth  rates  of  the  second  mode  of  instability  (in  particular). 
At  an  early  stage  it  began  very  clear  that  this  is  indeed  readily  possible,  and 
a  useful  choice  of  cj)  in  this  respect  is  =  7t/2.  For  this  reason,  the  results 
shown  here  will  be  confined  to  this  particular  choice,  and  emphasis  will  be 
place  on  the  effect  of  varying  the  ratio  . 

Figure  1  shows  results  for  a  square  cavity  (H  =  1),  at  a  Reynolds  number 
Re  =  20  and  a  frequency  parameter  of  0  =  1.  It  is  evident  from  these  results 
that  a  increase  in  the  parameter  serves  to  increase  the  magnitude  of 
both  -R{A* }  and  of  A* },  and  is  thereby  likely  to  improve  (in  practical 
terms)  the  stability  characteristics  of  the  outer  flow.  This  trend  is  repeated 
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Figure  4:  Singular  case,  £1  —  1,  Re  =  20,  ([>  =  n/2  (rectangular  cavity, 
H  =  2) 


in  the  results  shown  in  Figure  2  for  which  all  parameters  remain  unchanged, 
except  that  the  frequency  parameter  =  5  in  this  case.  We  see  in  these 
results  a  further  encouraging  trend,  namely  an  increase  in  the  magnitude  of 
both  5R{A*}  and  of  9{A*}  (compared  with  the  =  1  results),  whilst  this  is 
observed  again  as  increases.  One  important  aspect  here,  with  respect  to 
the  improvement  with  an  increase  in  the  value  of  the  frequency  parameter 
is  that  in  linking  the  outer  flow  regime  to  the  microcavity  flow  regime, 
the  net  effect  of  scaling  will  be  to  increase  the  value  of  0  in  the  microcavity 
flow  regime.  Yet  another  likely  improvement  in  flow  characteristics  is  to 
be  observed  in  Figure  3;  in  this  case,  the  parameters  taken  are  a  repeat  of 
those  in  Figure  1,  except  that  the  Reynolds  number  Re  =  40.  In  this  case 
there  is  an  approximate  doubling  of  the  magnitude  of  A*  (compared  with 
the  Re  =  20  results) . 

Figure  4  show  the  results  corresponding  to  Figure  1  but  for  the  case 
of  a  cavity  of  aspect  ratio  two  (i.e.  H  =  2,  corresponding  to  a  deeper 
cavity).  What  is  particularly  remarkable  with  these  results  ids  that  the 
results  compare  extremely  closely  with  those  for  the  unity  aspect  ratio  case 
(i.e.  H  =  1)  as  depicted  in  Figure  1. 
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Vtop 


Figure  5:  Regularized  case,  0  =  1,  Re  =  20,  <fi  =  tt/2  (square  cavity) 


Results  -  model  (b) 

In  this  case,  because  the  x  velocity  component  on  the  top  boundary  at 
x  =  0  and  x  =  1  is  zero,  then  the  singularity  encountered  at  these  two 
locations  with  model  (a)  is  no  longer  an  issue  (indeed,  this  type  of  approach 
is  frequently  referred  to  as  the  regularized  problem) . 

Figure  5  shows  results  for  a  square  cavity  [H  =  1)  with  a  Reynolds 
number  Re  =  20  and  a  frequency  parameter  0  =  1  and  (as  in  all  previous 
results)  a  phase  shift  parameter  <j>  =  tt/2.  There  is  evidence  here  of  9{A*} 
becoming  positive  as  the  ratio  is  increased,  although  (again  in  practical 
terms)  the  net  result  of  increasing  this  parameter  is  likely  to  still  have  ben¬ 
eficial  consequences,  on  account  of  the  increase  in  the  magnitude  of  ^{A*}. 
Figure  6  corresponds  to  the  case  0  =  5  (all  other  parameters  unchanged 
from  Figure  6).  In  this  case  the  results  bear  a  very  close  resemblance  to 
those  presented  above  for  the  corresponding  singular  case  (Figure  2).  This 
in  itself  is  encouraging,  indicating  a  degree  of  universality,  independent  of 
our  flow  model.  This  observation  is  repeated  in  Figure  7,  which  repeats  the 
calculation  of  Figure  2,  but  with  the  higher  Reynolds  number  Re  =  40;  this 
indicates  a  good  deal  of  similarity  with  the  corresponding  singular  results 
presented  in  Figure  3. 
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3  Conclusions  and  future  directions 


This  study  has  clearly  indicated  the  feasibility  of  significantly  stabilizing 
a  hypersonic  boundary  layer  through  the  use  of  microcavities.  This  effect 
is  again  seen  to  be  robust  -  achievable  over  a  broad  regime  of  parameters 
(in  particular  from  the  results  for  the  porosity  parameters  A  (A*)  in  our 
analysis),  which  indicates  the  potential  for  a  real  practical  tool  in  laminar 
flow  control. 

One  particularly  encouraging  observation  in  our  results  is  that  there  is  a 
good  deal  of  similarity  between  the  results  for  the  singular  cavity  (model  (a) ) 
and  those  for  the  regularized  cavity  (b)).  The  conclusion  from  this  is  then 
that  there  is  some  independence  of  results  on  modeling,  adding  credence 
to  the  results.  A  further  positive  feature  of  the  study  is  the  indication 
that  the  results  are  quite  independent  of  the  dimensions  of  the  microcavity 
(compare  Figure  1  for  a  unit  aspect  ratio  cavity  with  Figure  4  for  an  aspect 
ratio  of  two).  This  again  points  to  the  desired  effects  of  the  porous  surface 
being  achievable  without  the  need  for  delicate  fine  tuning  of  the  wall  surface 
characteristics. 

The  results  presented  in  this  report  suggest  that  an  improvement  may 
be  achieved  in  outer  flow  stability  characteristics  with  the  following: 

•  An  increase  the  in  frequency  parameter  fl;  as  noted  already,  the  ‘match¬ 
ing’  between  the  outer  flow  and  the  microcavity  flow  (taking  into  ac¬ 
count  the  different  scalings  in  these  two  zones)  will  certainly  have  this 
effect. 

•  An  increase  in  the  microcavity  Reynolds  number  ( Re) .  i.e.  with  an 
increase  in  the  flow  speed  within  the  cavity  (although  this  could  obvi¬ 
ously  have  a  corresponding  adverse  effect  on  the  boundary-layer  flow) . 

•  An  increase  in  the  ratio  parameter  indeed  this  is  not  incompatible 
with  the  comment  above  with  regard  to  the  Reynolds  number. 

On  account  of  the  significant  number  of  physical  parameters  occurring 
in  the  full  problem,  only  limited  (albeit  potentially  important,  interesting 
and  useful)  regions  of  parameter  space  have  been  explored;  obviously  a  more 
detailed  investigation  could  reveal  other  important  effects,  and  a  wider  study 
could  used  in  a  more  detailed  optimization  process  in  order  to  stabilise,  as 
far  as  possible,  the  boundary-layer  flow. 

The  modeling  in  this  work  has  relied  on  a  number  of  fundamental  as¬ 
sumptions,  and  the  potential  exists  for  relaxing  a  number  of  these  in  turn. 
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One  direction  of  significant  interest  would  be  to  investigate  the  role  of  Stokes 
layers  that  must  inevitably  occur  in  the  stability  regime  of  interest.  The 
stability  analysis  undertaken  in  this  project  has  been  based  on  inviscid  sta¬ 
bility  analyses,  which  is  entirely  correct  and  appropriate  in  the  hypersonic 
flow  regime.  However  this  analysis  also  requires  the  inclusion  of  thin  Stokes 
layers  very  close  to  the  body  surface  which  (in  mathematical  terms)  are  re¬ 
sponsible  for  ensuring  the  no-slip  condition  is  satisfied  on  the  body  surface. 
Inspection  of  the  parameters  in  the  experiments  of  Rasheed  et  al  [2]  reveal 
that  the  cavity  dimensions  are  comparable  to  the  thickness  of  the  Stokes  lay¬ 
ers,  and  there  therefore  exists  the  potential  for  significant  flow  interaction 
between  the  microcavities  and  these  Stokes  layers. 

A  further  effect  that  would  be  of  great  interest  would  be  the  cavity/cavity 
interaction.  In  the  analysis  performed  to  date,  single,  isolated  microcavities 
have  been  the  main  focus  of  attention  (although  the  the  associated  work 
of  Theofilis,  some  attention  has  been  directed  at  the  interaction  between 
just  two  cavities).  This  effect  is  certainly  an  area  which  could  shed  further 
significant  light  on  the  nature  of  the  overall  stabilization  process. 
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Recent  theoretical  and  experimental  investigations  have  suggested  that  the  mode  II 
instability  in  hypersonic  boundary  layers  can  be  considerably  stabilized  by  coating  the 
wall  surface  with  a  porous  material,  effectively  introducing  an  array  of  microcavities  on 
this  surface.  This  paper  is  concerned  with  two  aspects  of  this  problem.  Firstly  the  flow 
inside  a  solitary  microcavity  is  examined.  Order-of- magnitude  arguments  lead  to  the 
conclusion  that  the  flow  inside  each  microcavity  is  likely  to  be  not  only  incompressible, 
but  also  has  a  (local)  Reynolds  number  which  is  modestly  low;  this  is  due  to  the  small  scale 
of  the  microcavities.  Both  the  basic  flow  within  a  cavity  and  its  stability  are  calculated 
in  the  light  of  this  conclusion.  Secondly  the  flow  outside  the  microcavities  is  examined. 
A  flat-plate  hypersonic  boundary  layer  is  assumed  for  the  basic  state,  and  an  (inviscid) 
stability  analysis  (which  is  certainly  appropriate  for  mode  II  instabilities  in  the  Mach 
number  regime  of  interest)  is  performed.  The  porosity  of  the  wall  surface  in  modeled  in 
the  stability  analysis  by  permeability.  This  indicates  that  it  is  extremely  easy  to  halve 
disturbance  temporal  growth  rates  over  a  broad  range  of  parameter  values,  thus  indicating 
a  potentially  robust  and  hence  useful  means  of  flow  stabilization. 
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It  has  been  recently  shown  by  Fedorov  et  al.2 
theoretically  and  Rasheed  et  al.6  experimentally 
that  mode  II  instability,  which  prevails  in  hypersonic 
boundary-layer  flow  under  a  variety  of  environmental 
conditions,  can  be  very  effectively  controlled  passively 
by  coating  the  surface  on  which  the  instability  develops 
by  a  porous  material.  The  key  feature  of  this  technol¬ 
ogy  is  the  large  disparity  of  scales  between  a  typical 
mode  II  wavelength  and  the  diameter  d  of  each  of  the 
microscopic  cavities  on  the  coating,  schematically  de¬ 
picted  in  Figure  1.  One  consequence  is  that  the  effect 
of  an  array  of  microcavities  on  the  flow  can  be  mod¬ 
eled  by  replacing  the  impermeability  and  temperature 
boundary  conditions  at  the  wall  by  alternative  con¬ 
ditions  which  these  previous  authors  derived  partly 
under  the  assumption  of  infinitely  deep  cavities.  While 
this  approach  provides  a  satisfactory  first  estimate  of 
the  effect  of  the  coating  on  the  flow,2  potential  for 
improvement  of  the  performance  of  the  coating  exists 
through  improved  knowledge  of  the  flowfield  in  the 
neighborhood  of  each  individual  microcavity,  which 
would  thereby  permit  optimizing  the  spacing  S  be¬ 
tween  the  microcavities,  together  with  the  depth  D  of 
these  cavities,  both  of  which  are  currently  chosen  in  an 
arbitrary  manner.  Experiments  which  will  partly  fo¬ 
cus  on  this  issue  are  underway  in  Novosibirsk,  although 
their  expensive  nature  suggests  that  prior  guidance  of¬ 
fered  by  analytical/numerical  means  would  contribute 
to  focusing  the  efforts  on  selected  promising  configu¬ 
rations.  This  work  is  aimed  at  providing  the  first  esti¬ 
mate  of  the  necessary  information  by  numerical  means, 
together  with  an  improvement  in  our  scientific  under¬ 
standing  of  the  processes  that  occur  inside  each  cavity. 
We  consider  two  distinct  aspects  of  the  problem.  The 
first  aspect  is  concerned  with  developing  a  thorough 
understanding  of  the  flow  regime  inside  the  microcav¬ 
ities  themselves,  in  order  that  accurate  flow  modeling 
and  stability  analysis  can  be  performed.  The  second 
aspect  is  a  detailed  investigation  of  the  flow  instabil¬ 
ities  outside  of  the  microcavities,  i.e.  the  flow  in  the 
outer  boundary  layer.  This  study  leads  to  some  very 
interesting  conclusions  regarding  the  mimimization  of 
growth  rates,  with  obvious  application  to  the  practical 
problem  which  is  the  subject  of  this  project.  The  ul¬ 
timate  goal  will  be  to  link  the  two  flow  regimes  (that 
inside  the  microcavities  with  that  outside)  in  order  to 
predict  the  best  conditions  for  flow  stabilizations.  Our 
results  suggest  that  significant  flow  stabilization  can 
be  achieved,  over  a  broad  range  of  parameter  space, 
indicating  a  significant  and  feasible  practical  tool. 

2.  Steady  flow  inside  a  microcavity 


As  a  first  step  in  developing  a  detailed  model  for  use 
in  conjunction  with  the  (stability)  analysis  of  the  type 
described  below,  and  also  employed  by  Fedorov  et  al,2 


Fig.  1  Outline  of  the  problem  and  definition  of 
parameters 


i.e.  for  use  with  the  outer  flow  problem,  we  chose  to 
consider  the  flow  inside  a  single  two-dimensional  (mi¬ 
cro)  cavity.  Recent  work  by  Theofilis8  and  Theofilis, 
Duck  &  Owen9  has  considered  the  stability  of  the  lid- 
driven  flow  inside  a  rectangular  cavity  which  leads  to 
a  two-dimensional  eigenvalue  problem.  In  the  model 
under  consideration  here,  we  again  consider  the  basic 
flow  to  be  that  inside  a  rectangular,  two-dimensional 
cavity,  but  with  modified  boundary  conditions  at  the 
open  end  of  the  microcavity. 

The  theoretical  aspects  of  the  problem  decompose 
naturally  into  two  components:  the  basic  flow  within 
the  cavity  and  the  imposed,  small  amplitude,  distur¬ 
bance  field.  The  former  is  calculated  first,  and  then 
the  details  of  this  field  are  input  into  the  determina¬ 
tion  of  the  latter.  When  considering  the  base  flow,  a 
reasonable  starting  point  (given  the  spatial  dimension 
parameters  discussed  below)  is  to  impose  impermeabil¬ 
ity  together  with  uniform  shear  condition  on  the  cavity 
roof;  this  implies  that  to  leading  order  the  outer  flow  is 
unaffected  by  the  presence  of  the  microcavities.  On  the 
other  three  cavity  walls  the  usual  no-slip  and  imper¬ 
meability  conditions  are  obviously  appropriate.  These 
latter  conditions  are  clearly  also  entirely  correct  for 
the  flow  perturbations,  whilst  the  condition  appropri¬ 
ate  on  the  roof  of  the  cavity  for  the  disturbance  field  is 
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Fig.  2  Schematic  representation  of  the  minimal 
calculation  domain  encompassing  two  mi¬ 
crocavities 

less  clear.  Initially  zero  slip  and  impermeability  would 
seem  a  good  starting  point,  although  at  the  next  stage 
in  the  investigation,  given  the  results  to  be  described 
below,  the  impermeability  constraint  on  the  distur¬ 
bance  field  may  well  be  relaxed. 

In  the  light  of  the  remarks  above,  therefore,  just 
a  single  cavity  is  considered  in  isolation  -  given  the 
small  scale  of  the  cavity  width  (60  ±  4/<  m)  compared 
with  the  separation  between  cavities  (100  ±  7/i  m)  in 
the  experiments  of  Rasheed  at  al.6  this  would  seen 
to  be  a  reasonable  first  assumption;  in  turn,  both 
these  dimensions  are  considerably  smaller  than  the 
local  boundary-layer  thickness  in  the  outer  flow.  A 
further  underlying  simplification  is  that  the  basic  flow 
field  may  be  taken  as  being  two  dimensional,  a  pro¬ 
cedure  followed  in  the  theoretical  work  of  Fedorov  et 
al.,2  although  the  present  study  of  this  aspect  does 
allow  the  consideration  of  three-dimensional  imposed 
disturbances,  with  a  harmonic  decomposition  in  the 
spanwise,  i.e.  homogeneous  direction.  The  key  physi¬ 
cal  parameters  are  ,  the  freestream  Mach  number, 
the  thickness  of  the  incoming  boundary  layer  <5  (which 
is  uniquely  related  to  the  free-stream  velocity  U0 0)  one 
length  (e.g.  S  defined  in  Figure  2)  and  the  length  ra¬ 
tios  d/D  and  d/S,  as  well  as  fluid  viscosity  fx.  In  brief, 
therefore,  the  basic  flow  model  is  that  of  flow  inside 
a  rectangular  cavity,  with  no-slip/impermeable  sides 
and  base,  with  an  impermeable  top  with  an  imposed 
shear  stress.  With  this  model  it  will  be  possible  to 
refine  and  relax  a  number  of  the  above  assumptions. 

Of  paramount  importance  to  the  theoretical  work  is 
that  the  size  of  the  microcavities  is  such  that  we  may 
expect  that  the  velocity  inside  the  cavity  to  be  small, 
in  spite  of  the  hypersonic  nature  of  the  external  flow. 
This  has  two,  key  implications:  firstly  the  Reynolds 
numbers  based  on  the  microcavity  geometry  are  quite 
small  (indeed  the  scale  of  the  cavity  width  reinforces 
this) ;  secondly,  we  certainly  expect  that  the  flow  inside 
the  cavity  will  be  effectively  incompressible  (although 
temperature  variations  could  be  progressively  intro¬ 
duced  into  the  proposed  model  at  a  later  stage).  Let 
us  be  more  specific,  since  these  are  crucial  points.  In 
particular,  suppose  we  take  the  wall  shear  stress  of 
the  outer  flow  to  be  tw  =  /xw  ^ \w ,  then  we  may  expect 


freestream  and  wall  values  respectively.  If  the  cav¬ 
ity  (basic)  flow  is  triggered  in  non-dimensional  terms 
with  a  shear  stress  of  unity,  then  matching  this  quan¬ 
tity  inside  the  cavity  with  the  external  (wall)  value 
leads  us  immediately  to  the  conclusion  that  a  typi¬ 
cal  velocity  inside  the  cavity  is  U cavity  =  flToo,  i.e. 
|  times  the  freestream  value.  In  the  results  of  Fe¬ 
dorov  et  al.,2  this  translates  to  a  cavity /freestream 
ratio  of  0(0.06);  this  in  fact  assumes  that  the  typi¬ 
cal  cavity  basic  flow  velocity  is  0(1)  with  the  imposed 
unity  roof  shear  stress,  but  the  basic  flow  calculations 
performed  indicate  that  this  turns  out  to  be  consid¬ 
erably  smaller,  further  validating  our  argument.  The 
repercussions  for  the  Reynolds  number  are  even  more 
startling,  and  lead  to  the  conclusion  that  the  typical 
cavity  Reynolds  number  (based  on  width)  is  O(ys-) 
that  of  the  freestream.  In  practical  terms,  taking  the 
Fedorov  et  al.2  parameters,  this  translates  to  a  factor 
0(0.0036)  times  that  of  the  external  flow  or  a  Reynolds 
number  of  0(14)  taking  their  case  of  an  external  value 
of  4000.  This  strongly  suggests  that  the  Fedorov  et 
al.2  approach  of  using  inviscid  theory  to  predict  the 
behavior  of  the  flow  inside  each  cavity  may  not  be  suf¬ 
ficient  -  our  order-of-magnitude  arguments  suggest  a 
very  viscous  cavity  flow. 

A  number  of  simplifications  can  be  employed  in  or¬ 
der  to  gain  first  insights  into  the  configuration  at  hand. 
Those  pertaining  to  the  external  boundary-layer  flow 
are  discussed  below;  here  we  concentrate  on  the  sim¬ 
plifications  of  the  flowfield  inside  the  microcavities. 
These  amount  to  neglecting  the  three-dimensionality 
of  the  geometry  for  the  purposes  of  the  present  effort 
as  far  as  the  basic  state  is  concerned.  Regarding  global 
instabilities,  a  Fourier  decomposition  of  the  spanwise 
spatial  direction  is  considered.  As  has  been  mentioned, 
it  has  been  argued  that  flow  inside  the  microcavity 
is  essentially  incompressible  driven  by  constant  shear, 
which  leaves  two  parameters  to  be  considered  besides 
Reynolds  number,  namely  the  length  over  the  depth  of 
each  microcavity  and  the  spacing  between  two  neigh¬ 
bouring  microcavities.  In  the  absence  of  prior  infor¬ 
mation  on  the  flow  at  hand,  here  we  focus  on  the  basic 
state  inside  a  single  microcavity  and  its  global  insta¬ 
bility.  The  latter  is  contrasted  with  that  in  the  classic 
lid-driven  cavity.8  Along  the  N-boundary  of  the  inte¬ 
gration  domain  uniform  shear  has  been  imposed,  while 
no-slip  conditions  have  been  used  on  the  solid  walls  to 
close  the  problem. 

The  basic  flow  is  calculated  by  numerical  solution  of 
the  vorticity  transport  equation  alongside  the  relation 


between  streamfunction  and  vorticity, 

^  +  < 

f  dip  d( 

=  o, 

(1) 

dt 

{  dy  dx 

dx  dy  J  S 

K^  +  WeC  =  °’ 
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dv/dx  is  the  vorticity  of  the  basic  flow  and  ip  is 
its  streamfunction,  for  which  u  =  dip/dy  and  v  = 
—dipjdx  holds.  An  accurate  and  efficient  direct  nu¬ 
merical  simulation  algorithm  for  the  solution  of  (1)  - 
(2)  is  presented  in  the  next  subsection. 

The  linear  disturbance  equations  governing  the 
global  linear  eigenvalue  problem  are  obtained  by  sub¬ 
stituting  a  decomposition  into  two-dimensional  0(1) 
basic  quantities  q  =  (u,  v,  0,p)T  and  three-dimensional 
0(e)  disturbance  quantities  q  =  q (x,y)el  {Pz~m} ^ 
with  q  =  ( u,v,w,p)rT  into  the  equations  of  motion, 
subtracting  out  the  basic  flow  terms  and  neglecting 
terms  at  0(e2).  In  the  present  temporal  framework 
we  take  /3  to  be  a  real  wavenumber  parameter  de¬ 
scribing  an  eigenmode  in  the  z— direction  and  solve 
for  the  complex  eigenvalue  fl.  Qr  =  the  real 

part  of  the  eigenvalue  is  related  with  the  frequency  of 
the  global  eigenmode  while  the  imaginary  part  is  its 
growth/damping  rate;  a  positive  value  of  fij  =  3{f2} 
indicates  exponential  growth  of  the  instability  mode 
q  in  time  t  while  fl;  <  0  denotes  decay  of  q  in  time. 
When  the  wavenumber  vector  is  perpendicular  to  the 
plane  on  which  the  basic  flow  (u,v,0,p)  develops,  as 
the  case  is  in  the  present  basic  flow,  it  is  possible  to 
simplify  the  two-dimensional  eigenvalue  problem  by 
rewriting  it  as  one  with  real  coefficients,  which  requires 
half  the  amount  of  storage  compared  with  the  stan¬ 
dard  form.9  The  absence  of  a  basic  flow  ^-velocity 
component  in  the  linear  operator  in  conjunction  of  the 
redefinitions 

Cl  =  i  fi,  w  =  i  w  (3) 

result  in  the  following  generalized  real  nonsymmetric 
partial  derivative  EVP  after  linearisation  and  subtrac¬ 
tion  of  the  basic-flow  related  terms, 


[M  -  (Dxu)\  u  -  (Vyu)v  -  Vxp 

=  fi  u. 

(4) 

~(Vxv)u  +  [M-  (T>yv)\  v  -  T>yp 

<sT 

II 

(5) 

+Mw  -  Bp 

=  Cl  w, 

(6) 

Vxii  +  Vvv  —  Bw 

=  0, 

(7) 

where 

M  =  (l/Re)  (Vl+V2y+p2)-uVx 

~VVy. 

(8) 

The  ability  to  solve  a  real  eigenvalue  problem  is  essen¬ 
tial  as  the  Reynolds  numbers  increases,  Indeed,  this 
point  has  been  clearly  manifested  in  the  literature  in 
the  difficulties  encountered  by  the  investigators  who 
tackled  the  problem  of  linear  instability  in  the  classic 
lid-driven  cavity  to  produce  consistent  results. 
Numerical  methods 

The  principles  of  the  algorithm  discussed  in 
what  follows  can  be  applied  to  recover  three  two- 
dimensional  velocity  components  or  indeed  a  three- 
dimensional  flowfield;4  however,  for  simplicity  we  con¬ 
fine  the  present  discussion  to  solutions  of  the  system 


main  advantage  of  the  velocity-vorticity  formulation  is 
that  the  continuity  equation  is  exactly  satisfied.  How¬ 
ever,  the  problem  of  imposition  of  boundary  conditions 
within  the  framework  of  an  overall  efficient  numerical 
solution  algorithm  remains.  This  is  compounded  by 
the  fact  that  the  number  of  points  discretising  the  two 
spatial  directions  in  the  subsequent  analysis  cannot  be 
increased  at  will;  while  interpolation  of  a  basic  flow  so¬ 
lution  obtained  on  a  very  large  number  of  points  onto 
a  modest  EVP  grid  is  one  possible  option,  it  is  more 
elegant  to  avoid  the  interpolation  procedure  altogether 
and  seek  an  accurate  basic  flow  solution  on  the  same 
small  number  of  discretisation  points  on  which  the  sub¬ 
sequent  global  instability  analysis  is  to  be  performed. 
This  appears  tailor-made  for  a  spectral  numerical  so¬ 
lution  approach;1  in  what  follows  we  discuss  a  different 
solution  based  on  an  efficient  real-space  eigenvalue  de¬ 
composition  (EVD)  algorithm. 

Chebyshev  polynomials  have  almost  exclusively 
been  used  in  the  past  in  the  context  of  spectral  simu¬ 
lations  of  the  time-accurate  Navier-Stokes  and  con¬ 
tinuity  equations,  mainly  due  to  the  availability  of 
fast  transform  algorithms  necessary  for  efficient  time- 
integration.  However,  for  the  present  problems  we 
have  not  restricted  ourselves  to  this  class  of  orthog¬ 
onal  polynomials.  Considerable  freedom  exists  in  the 
choice  of  the  expansion  functions  and  the  associated 
collocation  grids  by  using  Jacobi  polynomials  P^’O  for 
the  discretization  of  both  spatial  directions;  of  course, 
q  =  r  =  —0.5  may  be  related  to  the  Chebyshev-  while 
q  =  r  =  0  are  the  Legendre  polynomials.  Collocation 
derivative  matrices  for  both  Jacobi  Gauss-Lobatto  and 
equidistant  grids  can  be  constructed  from  first  princi¬ 
ples;  if  ( Xj ,  j  =  0,  ■  ■  ■ ,  n)  is  the  collocation  grid  chosen, 
the  entries  dij  of  the  (0  :  n)  x  (0  :  n)  first-order  deriva¬ 
tive  collocation  matrix  V,  derived  analytically  from 
the  interpolating  polynomial,  are 

n 

n  (xi  -  xk) 

dij  =  - k=0  - - ,  i.j.k  =  0,  -  ■  -,n. 

fa  -Xj)  n  (xj  - xk) 

k= 0 

(9) 

da  =  — - - - ,  i,k  =  0,---,n,  i  ±  k.  (10) 

E  (xi  ~  xk) 

k=  0 

These  formulae  result  in  the  well-known  ones  if  the  an¬ 
alytically  known  Chebyshev  Gauss-Lobatto  grid  (xj  = 
cos  jn/n,  j  =  0,  ■  ■  ■  ,n)  is  used.1  Values  of  order  m 
derivatives  on  the  collocation  grid  Xj  are  obtained  by 
(' D)m .  As  far  as  temporal  discretization  of  (1)  is  con¬ 
cerned,  the  viscous  nature  of  the  problems  in  which 
we  are  interested  introduces  scales  which  dictate  an 
implicit  treatment  of  the  linear  term  in  this  equa¬ 
tion;  the  nonlinear  term  may  be  treated  explicitly. 
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not  resort  to  splitting  and  introduction  of  intermedi¬ 
ate  fields  but  rather  address  the  governing  equations 
directly  the  combination  of  Crank-Nicholson  (CN) 
with  second-order  Adams-Bashforth  (AB2)  or  Runge- 
Kutta  (RK)  schemes  has  been  extensively  used  for 
the  time-integration  of  the  viscous  and  the  convective 
terms,  respectively.1  However,  the  family  of  com¬ 
pact  schemes  proposed  by  Spalart,  Moser  and  Rogers7 
(SMR)  presents  more  accurate  and  more  stable  alter¬ 
natives  to  the  CN-AB2  algorithm  although  it  does  not 
require  additional  computational  effort  to  the  latter 
scheme.  The  SMR  algorithm  may  be  written  in  com¬ 
pact  form  as 

«2  W(g')}, 

(11) 

where  the  superscript  denotes  fractional  time-step, 
C(q)  and  J\f(q)  are,  respectively,  the  linear  and  non¬ 
linear  operators  in  the  problem  to  be  solved  and  At 
is  the  time-step.  The  rationale  behind  the  derivation 
as  well  as  sample  values  of  the  constants  k.  A,  K\  and 
K2  of  a  self-starting  SMR  algorithm  may  be  found  in 
Spalart  et  al.7  The  time-discretization  of  (1)  using 
(11)  delivers  the  following  problem  to  be  solved  for 
((.  ip)  at  each  fractional  time-step 


Mi  (  =R, 

(12) 

ill  ill 

M2  ip  +  C  =  o, 

(13) 

where  Mi  =  d2  jdx2  +d2/dy 2  -  Re/(XAt)  and  M2  = 
d2 /dx2  +  d2 / dy2,  subject  to  boundary  conditions  ap¬ 
propriate  to  the  problem  in  consideration.  R  comprises 
the  nonlinear  and  the  terms  arising  from  the  discretiza¬ 
tion  at  previous  fractional  time-steps, 


rithm  may  be  constructed  for  the  efficient  solution 
of  the  incompressible  two-dimensional  Navier-Stokes 
and  continuity  equations  in  streamfunction-vorticity 
formulation. 

Key  papers  on  EVD  algorithms  are  the  work  of 
Haidvogel  &  Zang3  and  that  of  Ku  et  al.4  The  first 
authors  discuss  the  solution  of  Poisson’s  eq 


III  ll  f  "  ///  .  II 

q  =q  +  AtijC(nq  +  Xq  )  +  KiXf(q  ) + 


R  = 


+ 

+ 


-("/A) 


d 2  d2 

ft? +  -*/<«*) 


(tiiRe/X)  Cy) 

(K2Re/X)  (ip'X  -  ipx' 4)  •  (14) 


The  accuracy  of  the  overall  procedure  clearly  de¬ 
pends  on  the  scheme  utilised  for  calculation  of  the 
spatial  derivatives.  The  spectral  discretisation  cho¬ 
sen  introduces  dense  matrices  and  can  only  become 
competitive  against  other  numerical  approaches  from 
the  point  of  view  of  efficiency  on  account  of  the  ex¬ 
istence  of  a  fast  algorithm  for  the  inversion  of  the 
implicit  operators  Mi  in  (12)  and  M2  in  (13).  While 
M2  is  time-independent,  the  first  implicit  operator 
is  a  function  of  a  CFL-controlled  maximum  permit¬ 
ted  time-step  At  and  needs  to  be  inverted  at  every 
time-step.  If  one  sacrifices  the  advantage  of  an  ad¬ 
justable  time-step  and  keeps  At  fixed  at  a  slightly 


On  account  of  the  homogeneous  boundary  conditions 
(15)  on  i  (13)  becomes 

m— 2  n— 2 

y  Mkiipu  +  y  Nijipkj  =  -  (ki  -  Mki'tpu 

i=  2  j=2 

-  Nin-itpkn-1  ~  Nnlpkl 
—  Mkm  —  llpm  —  ll-  (18) 

The  boundary  conditions  (16)  -  (17)  may  be  expressed 
using  the  discrete  analoga  X,  Y  and  Y2  of  the  colloca¬ 
tion  derivative  matrices  Vx,  Vy,  and  V2,  respectively, 
as  given  by  (9)  -  (10).  It  follows  that 


V’fci 


m— 2 


=  y  Sain 

i= 2 
n- 2 

=  Kk  +  y  ei jikj 
1=2 


'Ipm—ll 


'Ipkn—l 


m— 2 

^  ^  (19) 

i=2 


n- 2 

“1“  ^  ^  €n—lj'*ftkj  i 
3= 2 


diagonalised  using  their  eigenvalue  decomposition 
M  =  (M*)^*(M*)-1  and  iV  =  (29) 

to  yield 

^(M*)-1/W(1V*)+(M*)-1/«(1V*)  =  (M*)-VW(1V*). 

(30) 

As  a  consequence,  instead  of  having  to  solve  the 
(to  —  3)  x  (n  —  3)  system  of  simultaneous  equations 
(28)  one  solves  the  (to  —  3)  x  (n  —  3)  algebraic  equa¬ 
tions  /*  =  r*/(n*+v*)  for  /*  =  (AT)”1  /H(JV'),  given 
r*  =  ( M*)~1rki(N *).  The  cost  of  the  algorithm  is  a 
negligibly  small  fraction  of  the  cost  of  a  direct  algo¬ 
rithm  for  the  solution  of  the  Poisson  problem.  This 
is  documented  in  table  1  where  memory  and  run-time 
requirements  are  shown  for  solution  of  (1)  -  (2)  in  the 
lid-driven  cavity  problem. 

3.  Results 


(20) 

where  Rk  and  A*,  are  used  to  impose  the  boundary 
condition  on  the  roof  of  the  microcavity, 


_  _  -FkYon-i  -r  _  FkY0 1 

Kk  ~  v  v 2  _  V2  v  ’  —  A''  v2  —  v2  v 

10llnn-l  1  nl1  On-1  i01-iTOi_i  *  nl^  On— 1 

(2i) 

and  the  vectors  Sii,dm-u  and  eij,  en-ij  are  known 
functions  of  the  entries  of  X,  Y  and  Y2, 


$li 


<5 


m—li 


£i  j 


en—lj 


—  Xmm-iXoi 


XmiXoui—i 

XmlXoi 

X01Xmm-l 

X-miXom—i 

Yon-lY^j  - 

V2  V  • 

1  nn—11  Oj 

YoiYL-i  ~ 

Y^Yon-!  ’ 

YniYoj  ~ 

Y0iY% 

YoiY2^  - 

Y^Yon-1  ■ 

(22) 

(23) 

(24) 

(25) 


The  essence  of  the  EVD4  algorithm  is  to  diagonalise 
the  (0  :  to  —  4)2  matrix  M  and  the  (0  :  n  —  4)2  matrix 
N  whose  entries  are 


—  Mki  T  MfciSn  + 

k, i  =  0,  ■  ■  ■  ,to  —  4,  (26) 

=  Nij  +  Nneij  +  Nin- ien-ij, 

l, j  =  0,  to -4.  (27) 

This  is  in  contrast  to  the  EVD2  algorithm  of  Ku  et 
al 4  where  an  (0  :  m  —  2) 2  matrix  M  and  an  (0  :  n  — 
2) 2  matrix  N  arising  from  the  imposition  of  Dirichlet 
data  must  be  diagonalised.  The  Poisson  problem  (18) 
becomes 


Mki 

Nij 


Mfki  +  fkiN  =  —(ki  ~  NnKk  -  Nin-iXk  =  rki  (28) 


Basic  flow 

We  implemented  the  algorithm  discussed  in  the  pre¬ 
vious  section  which  was  modified  compared  with  that 
of  Theofilis,  Duck  &  Owen9  by  replacing  the  uniform 
tangential  velocity  of  the  cavity  roof  by  that  of  uni¬ 
form  wall  shear  stress,  as  detailed.  The  issue  of  the 
singularity  of  the  boundary  conditions  at  the  NE  and 
NW  corners  of  the  cavity  is  thus  absent  in  the  present 
calculations  and  the  time-accurate  spectral  colloca¬ 
tion  scheme  utilized  for  the  calculation  of  the  basic 
flow  demonstrates  exponential  convergence.  Unlike 
the  standard  lid-driven  cavity,  in  which  the  condition 
iy  =  1  is  imposed  at  the  N  boundary  and  defines  the 
effective  flow  Reynolds  number,  there  is  some  ambi¬ 
guity  in  the  definition  of  a  Reynolds  number  in  the 
present  case.  The  input  Reynolds  number  may  be  re¬ 
garded  as  Re  =  uv(y  =  D)d2 /v,  given  that  in  our 
approach  the  roof  shear  stress  is  constant  (or  ivv  =  1); 
here  v  is  the  kinematic  viscosity  of  the  fluid.  Two  addi¬ 
tional  Reynolds  numbers  useful  for  closer  comparisons 
with  the  standard  lid-driven  cavity  problem  may  also 
be  defined,  namely  an  integral  Reynolds  number 

1  fd 

Reint  =  ~  /  u(x,  y  =  D)  dx,  (31) 

v  Jx= 0 

as  well  as  that  based  on  the  maximum  streamwise  ve¬ 
locity  component  attained  at  y  =  D, 

7?emax  =  ^  max{u(x,y  =  D)}.  (32) 

While  v  is  a  fixed  parameter  in  the  basic  flow  calcu¬ 
lation  code,  u(x,  y  =  D)  is  unknown  a  priori  and  is 
determined  from  the  converged  basic-flow  field.  Tak¬ 
ing  d  =  D  =  1  basic  flow  results  have  been  obtained 
using  a  rectangular  grid  comprising  upwards  of  96  col¬ 
location  points  per  spatial  direction.  These  have  been 
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(Mb) 

time 

(sec) 
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(Mb) 

time 

(sec) 

size 

(Mb) 

time 

(sec) 

size 

(Mb) 

time 

(sec) 
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0.4 

4.4 

1.1 

3.7 

4.03 

0.03 

4.03 

0.1 

242 

0.5 

5.4 

3.6 

10.5 

5.03 

0.14 

6.03 

0.3 

322 

0.6 

6.6 

10.0 

39.7 

5.03 

0.25 

13.03 

1.1 

482 

0.8 

15.3 

46.9 

460.6 

6.03 

0.56 

48.03 

8.4 

642 

1.1 

31.8 

143.9 

5203.5 

6.03 

1.08 

140.03 

40.7 

962 

1.8 

143.3 

(*) 

(*) 

8.03 

2.48 

680.03 

417.4 

1282 

2.8 

523.1 

(*) 

(*) 

8.03 

4.41 

(*) 

(*) 

Table  1  Comparison  of  memory  and  run-time  requirements  for  a  single  solution  of  a  two-dimensional 
Poisson  equation  using  direct  inversion  and  the  EVD4  algorithm  on  one  processor  of  a  work¬ 
station  and  a  supercomputer.  Asterisks  denote  that  the  respective  problem  does  not  fit  in  the 
available  memory  on  the  workstation  or  that  it  cannot  be  solved  within  the  existing  batch  queue 
time-limit  on  the  supercomputer. 


compared  against  the  classic  lid-driven  cavity  flow,  the 
Reynolds  number  in  the  latter  case  taken  to  be  of 
0(Re int)  in  the  former.  As  an  example  a  comparison 
between  the  constant-shear  and  lid-driven  basic  states 
in  terms  of  the  streamfunction  ip,  the  velocity  compo¬ 
nents  u  =  ipy,v  =  —ipx  and  the  vorticity  (  is  shown 
in  Figure  3.  The  shear-driven  cavity  results  have  been 
generated  using  Re  =  6250,  which  yields  an  effective 
Reynolds  number  of  Re jnt  «  300.  Both  sets  of  the 
first  three  quantities  and  the  vorticity  of  the  shear- 
driven  flow  are  presented  by  ten  isolines  between  zero 
and  the  respective  maxima.  The  comparison  points 
to  a  general  qualitative  agreement  between  the  two 
model  flows.  That  the  input  Reynolds  number  (Re)  in 
the  shear-driven  case  is  substantially  higher  than  Re int 
shows  how  relatively  ineffective  the  constant  shear 
driving  mechanism  is.  In  anticipation  of  the  parameter 
range  necessary  to  produce  globally  unstable  flows,8,9 
constant-shear  basic  states  were  obtained  at  substan¬ 
tially  larger  values  of  the  unit  Reynolds  number,  as 
shown  in  Figure  4.  A  monotonic  increase  of  both  i?eint 
and  I?emax  as  l/v  increases  can  be  seen.  The  extent 
to  which  the  quantitative  differences  between  the  ba¬ 
sic  states  modify  the  stability  characteristics  of  the 
lid-driven  cavity  flow  is  examined  next. 

Global  linear  instability 

In  order  to  address  the  two-dimensional  partial 
derivative  eigenvalue  problem  appropriate  boundary 
conditions  must  be  specified.  These  are  the  straight¬ 
forward  viscous  boundary  conditions  on  the  distur¬ 
bance  velocity  components  on  the  three  solid  walls 
while  the  issue  of  boundary  conditions  on  the  N  bound¬ 
ary  of  the  shear-driven  cavity  is  expected  to  be  the  key 
to  link  the  present  configuration  with  that  prevailing 


of  prior  information,  calculations  reported  here  have 
been  performed  by  imposing  homogeneous  Dirichlet 
or  extrapolation  of  the  disturbance  velocity  compo¬ 
nents  from  the  cavity  interior;  appropriate  compat¬ 
ibility  conditions  for  the  disturbance  pressure  have 
been  used  to  close  the  system  to  be  solved  in  botli 
cases.  In  the  framework  of  a  temporal  global  instabil¬ 
ity  analysis  the  only  parameters  of  the  problem  are  the 
flow  Reynolds  number  Re  and  the  spanwise  periodicity 
length  Lz,  associated  with  a  wavenumber  B  =  27t/Z/z. 
while  the  complex  temporal  eigenvalue  fl  and  the  two- 
dimensional  amplitude  functions  {u,v,w,p)T  of  the 
disturbance  quantities  are  determined  by  the  numeri¬ 
cal  solution  of  the  partial-derivative  eigenvalue  prob¬ 
lem. 

At  Re ide  =  Reint  «  300  the  lid-driven  cavity  flow 
is  known  to  be  stable  to  both  two-  and  all  three- 
dimensional  perturbations.  Comparisons  between  the 
two  model  flows  are  only  possible  in  the  case  of  homo¬ 
geneous  Dirichlet  boundary  conditions  being  imposed 
at  the  N  boundary  in  the  shear-driven  case.  The 
modification  of  the  (least  damped)  leading  eigenval¬ 
ues  of  the  eigenspectrum  at  /?  =  7.5  (approximately 
corresponding  to  least-damped  conditions  for  mode 
T28,9)  is  shown  in  Figures  5  and  6.  Owing  to  the 
disparity  between  i?eidc  and  Re  the  shear-driven  flow 
experiences  much  weaker  damping  compared  with  its 
lid-driven  counterpart.  Noteworthy  is  also  the  differ¬ 
ence  between  the  two  flows  in  terms  of  the  prevalence 
of  stationary  and  low-frequency  disturbances  (in  the 
neighbourhood  of  =  0)  in  the  shear-driven  cavity 
case  at  Re  =  6250.  At  Re  =  25  x  104,  on  the  other 
hand,  the  flow  is  marginally  unstable  at  /?  =  7.5.  The 
amplitude  function  of  the  spanwise  disturbance  veloc¬ 
ity  component  w  of  the  unstable  global  eigenmode  is 
shown  in  Figure  7.  The  complexity  of  the  (periodic  in 
z)  flowfield,  setup  by  global  instability  inside  the  mi¬ 
crocavity,  can  be  seen  in  this  figure,  while  the  analogies 
of  the  structure  of  the  global  eigenmode  with  that  in 
the  lid-driven  cavity  can  also  be  appreciated.8,9 
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Fig.  4  Dependence  of  R<:\n\  and  Remax  on  1  / v  in  a 
square  shear- driven  cavity. 
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4.  The  flow  outside  the  microcavities 

Turning  our  attention  to  the  outer  boundary  layer  flow 
next,  the  basic  state  in  this  region  has  been  calcu¬ 
lated  using  the  standard  compressible  boundary-layer 
equations,  together  with  the  associated  energy  equa¬ 
tion,  assuming  the  similarity  form  appropriate  for  a 
boundary  layer  on  a  flat  plate  with  zero  streamwise 
pressure  gradient.  The  Chapman  viscosity  law  has 
been  utilised;  although  this  could  easily  be  substituted 
for  a  more  sophisticated  model,  our  study  indicates 
than  any  resulting  differences,  within  the  Mach  num¬ 
ber  range  considered,  are  minimal.  Note  also  that 
the  effects  of  shock  waves  have  been  ignored  (indeed, 
this  approximation  was  also  employed  by  Fedorov  et 


al.2).  We  note  that  the  shock  waves  in  the  Mach  num¬ 
ber  range  of  immediate  interest  are  expected  to  lie 
well  outside  of  the  boundary  layers,  and  that  the  flow 
disturbances  of  most  interest  are  anticipated  to  be  con¬ 
fined  to  well  inside  the  boundary  layers.  Results  for  the 
base  state  are  not  presented  here  -  these  are  well  doc¬ 
umented,  and  our  results  are  entirely  consistent  with 
such  previous  work. 

It  is  well  known  (Mack5)  that  an  inviscid  stabil¬ 
ity  analysis  of  high-supersonic/hypersonic  boundary 
layers  iba  tfe .  id  bbla  nip  P  t  a  |>  fe|» 
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Fig.  6  Global  eigenspectra  in  the  lid-driven  cav¬ 
ity. 


tainly  interested  in  here.  We  note  that  Fedorov  et 
al.2  did  perform  external  stability  calculations,  using 
the  full  viscous  disturbance  equations,  but  the  pro¬ 
cedure  adopted  here,  of  using  inviscid  stability  anal¬ 
ysis  permits  many  such  calculations  to  be  performed 
rapidly  (on  a  PC /laptop  computer) ,  and  thereby  per¬ 
mits  a  thorough  investigation  of  parameter  space,  and 
a  subsequent  optimization  process  to  be  performed  ef¬ 
ficiently. 

Specifically,  taking  the  x  coordinate  in  the  stream- 
wise  direction,  Y  the  normal  boundary-layer  coor¬ 
dinate  (scaled  with  Re~ i  in  the  usual  boundary- 
layer  manner) ,  then  the  corresponding  base  flow  (non- 
nondimensionalized  with  respect  to  the  freestream  ve¬ 
locity  U0 o,  which  is  also  used  in  our  definition  of  the 
Reynolds  number  Re),  and  temperature  field,  takes 
the  form 

(u,T)  =  ((Uo(Y),0),To(Y))  +0(Re~v).  (33) 

Notice  here  the  base  flow  is  taken  to  be  inde¬ 
pendent  of  the  streamwise  coordinate  x  -  in  reality 
such  a  locally  parallel  flow  analysis  is  perfectly  cor¬ 
rect/justifiable  if  the  (following)  stability  analysis  is 
based  on  disturbances  comparable  to  the  boundary- 
layer  thickness  -  this  is  certainly  the  case  here. 

We  therefore  decompose  the  perturbation  flowfield 
in  the  following  (normal-mode)  manner 


(u,v,P,T,p)  =  (u,v,p,f,p)eiR^a^~ct\ 


(34) 


where  u,  v  and  T,  take  on  obvious  meanings,  p 
and  p  represent  the  perturbation  pressure  and  den¬ 
sity,  respectively,  and  a  denotes  the  wavenumber  of 
the  disturbance,  c  its  wavespeed;  the  amplitude  func¬ 
tions  in  the  following  local  analysis  are  one-dimensional 
functions  of  the  wall-normal  coordinate  Y.  by  con¬ 
trast  to  the  global  analysis  of  the  previous  section, 
where  the  amplitude  functions  were  defined  on  the 
two-dimensional  plane  of  the  microcavity.  Note  also 
that  we  are  assuming  a  two-dimensional  disturbance 
field  -  as  shown  by  Mack  this  is  the  most  ’dangerous’ 
mode  in  this  regime;  however  it  is  a  simple  matter  to 
extend  the  analysis  into  the  three-dimensional  distur¬ 
bance  regime. 

Substituting  the  disturbance  flowfield  into  the  in¬ 
viscid,  compressible  equations  of  motion,  discarding 
quadratic  terms  leads  to  the  following  system 


Fig.  7  Four  isosurfaces  of  the  spanwise  distur¬ 
bance  velocity  inside  a  square  microcav¬ 
ity  at  (Re  =  25  x  104,  /3  =  7.5).  Levels  are 
equidistributed  between  the  minimum  and 
maximum  values. 
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P  =  J-OP  + 


To' 


After  some  algebra  it  is  possible  to  eliminate  all  but 
two  of  the  dependent  variables,  and  our  choice  was  to 
work  with  the  following  system 


Vy  ~ 


U0Yv 

U0-c 


ip  \T0  -  Ml(U0  -  c)2' 

1M^  [  U0-c  J  ’ 

ia2(Uo—c)v  pY 

T0  =~Wl' 


(40) 

(41) 


In  the  usual,  standard  case  of  zero  permeability  on 
the  wall  the  appropriate  the  boundary  conditions  are 
that  v(Y  =  0)  =  0,  together  with  (generally)  the  dis¬ 
turbance  field  decays  as  Y  — >  oo  *.  However  it  is  our 
assertion  (which  agrees  with  Fedorov  et  al.2)  that  the 
effect  of  the  microcavities  is  primarily  to  modify  this 
surface  condition.  We  choose  therefore  to  impose  the 
following  modified  condition,  which  replaces  imperme¬ 
ability  with  a  relationship  between  normal  velocity  and 
pressure  on  the  wall 


v  —  Ap  =  0. 


(42) 


This  is  a  reasonably  general  relationship,  consistent 
with  a  well-posed  solution  of  the  above  system.  A 
numerical  code  was  developed  to  solve  the  temporal 
eigenvalue  problem,  treating  a  as  a  real,  specified  pa¬ 
rameter  and  evaluating  the  corresponding  eigenvalue 
c.  The  method  was  based  on  a  fourth-order  solver, 
which  has  proved  fast  and  effective  on  problems  of  this 
type  in  the  past  (Mack  ).  In  the  case  of  neutral/stable 
disturbances,  it  is  necessary  to  perform  contour  in¬ 
dentations  (into  the  complex  Y  plane),  and  this  was 
an  integral  part  of  the  numerical  code.  A  very  thor¬ 
ough  and  detailed  study  was  then  made  over  a  range  of 
wavenumbers  (a)  and  of  surface  parameter  (A).  Our 
results  focus  on  the  temporal  growth  rate  0,  =  3{ac}. 
Note  there  is  no  reason  why  A  should  necessarily  be  a 
real  quantity,  and  so  we  therefore  conducted  a  detailed 
study  in  complex  A  space,  and  we  write  A  =  Ar  +  i A,. 
Further,  in  order  to  replicate  the  flow  conditions  stud¬ 
ied  by  Fedorov  et  al.2  and  Rasheed  et  al.6  all  results 
shown  here  correspond  to  a  freestream  Mach  number 

Moo  =  6. 

Figure  8  shows  results  for  A*  =  0,  for  selected  val¬ 
ues  of  Ar ;  the  case  Ar  =  0  therefore  corresponds  to  the 
standard  impermeable  results  (which  can  be  checked 
against  those  of  previous  studies,  for  example  Mack  ). 
In  this  standard  case,  the  first  and  second  modes  of  in¬ 
stability  are  quite  distinct,  and  in  line  with  the  usual 
result  that  it  is  the  second  mode  that  dominates  at 
these  high  Mach  numbers.  As  Ar  is  reduced,  however, 

’There  are  some  exceptional,  neutrally  stable  cases  where  the 
disturbance  is  bounded  as  t  — >  oo  -  we  do  not  present  details 
here,  but  we  do  note  that  our  numerical  scheme  is  well  able  to 
handle  these. 


a 

Fig.  8  A,:  =  0 


a 

Fig.  9  A,:  =  — 1 


initially  there  is  an  amalgamation  of  the  two  modes, 
and  also,  crucially  a  marked  reduction  in  maximum 
growth  rate.  In  this  case,  the  optimum  choice  of  Ar 
(i.e.  that  exhibiting  the  least  growth  rate)  is  close  to 
—4  and  is  shown  as  a  dashed  line  (this  is  the  con¬ 
vention  we  also  use  below  to  indicate  near  optimum 
conditions).  Thereafter,  i.e.  as  Ar  is  reduced  still  fur¬ 
ther,  an  increase  in  growth  rates  is  observed  at  higher 
values  of  a  (perhaps  an  increase  in  the  significance 
of  mode  III)  and  this  soon  dominates  and  leads  to  a 
worsening  of  the  stability  characteristics.  What  is  ap¬ 
parent,  however,  is  that  in  this  case,  i.e.  A*  =  0  it  is 
very  easy  to  get  an  approximate  halving  of  maximum 
growth  rate  over  a  very  broad  range  of  Ar  -  this  pa¬ 
rameter  does  not  require  accurate  or  sensitive  tuning 
to  achieve  this  result. 

Figure  9  shows  the  corresponding  results  when  A*  = 
—  1.  There  is  a  good  deal  of  similarity  displayed  with 
those  of  Figure  8,  although  in  this  case  it  is  the  choice 
Ar  «  —  2  that  exhibits  the  lowest  maximum  growth 
rates;  nonetheless  an  approximate  halving  of  maxi¬ 
mum  growth  rate  is  achieved  again  over  a  broad  range 
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Fig.  10  A  i  =  -2 


Fig.  12  A  i  =  -10 


a 

Fig.  11  A,  =  — 5 


we  may  expect  improved  modeling  of  the  flow  using 
the  more  realistic  (viscous)  equations. 

Some  detailed  stability  analysis  has  also  been  per¬ 
formed  on  a  boundary  layer  in  a  Mach  6  external  flow. 
Using  a  modified  surface  boundary  condition  for  this 
analysis  to  model  the  effects  of  the  microcavities,  it  has 
been  shown  vividly  that  it  is  straightforward  to  obtain 
a  halving  of  the  maximum  growth  rates  over  a  very 
wide  range  of  surface  conditions.  This  suggests  a  use¬ 
ful  (and  robust)  tool  for  transition  delay /prevention. 

Future  work  will  focus  on  use  of  the  boundary  con¬ 
dition  at  the  open  end  of  the  microcavity  as  the  key  to 
relate  flow  instability  inside  the  boundary  layer  with 
that  inside  the  microcavity  itself,  in  the  continuing 
quest  for  hypersonic  flow  control  through  control  of 
flow  instability. 


of  Ar. 

Results  for  =  —2  are  displayed  in  Figure  10.  In 
this  case  the  optimum  value  of  Ar  an  —2,  but  the  re¬ 
sults  are  otherwise  in  qualitative  agreement  with  those 
of  Figures  8  and  9.  Results  for  A,  =  —5  and  A,  =  —10 
are  shown  in  Figures  11  and  12,  respectively,  and  ex¬ 
hibit  optimum  stability  characteristics  at  Aj  «  —  7  and 
Aj  «  — 12  respectively.  Remarkably,  all  our  results  il¬ 
lustrate  how  it  is  a  simple  matter  to  approximately 
halve  maximum  (temporal)  growth  rates  over  a  wide 
range  of  choices  of  the  parameter  A.  This  strongly 
suggests  this  may  be  a  useful  and  quite  straightfor¬ 
ward  mechanism  for  outer  flow  stabilization. 

5.  Conclusions 

In  paper  we  have  established  the  nature  of  the  basic 
flow  inside  a  cavity,  driven  by  a  uniform  shear-stress 
on  one  face.  This  has  been  determined  to  be  primarily 
viscous  in  nature,  and  thereby  justifies  the  use  of  the 
(full)  viscous  and  incompressible  equations  in  this  re¬ 
gion,  as  employed  in  this  paper.  We  do  note  that  in  the 
earlier  work  of  Fedorov  et  al.2  inviscid  equations  were 
used  inside  the  microcavities  (a  detailed  comparison 
will  be  made  with  this  work  in  due  course)  although 
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